import numpy as np
from scipy.spatial.transform import Rotation

class sate_motion(object):
    def __init__(self, microthrust, sate_mass) -> None:
        '''
        计算卫星运动所需时间及速度增量
        microthrust: 微推推力，单位 N
        sate_mass: 卫星质量，单位 kg
        '''
        self.__microthrust = microthrust
        self.__sate_mass = sate_mass

        self.__max_acceleration = self.__microthrust / self.__sate_mass

    def linear_motion(self, distance):
        '''
        L = 0.5 * a * t^2
        "Bang-bang"机动模式
        '''
        delta_v = np.sqrt(self.__microthrust * distance / (9. * self.__sate_mass))
        delta_t = 10 * distance / (9*delta_v)

        return delta_t / 3600., delta_v*2

    def circling_motion(self, radius):
        #max_angvec = np.sqrt(self.__microthrust / (self.__sate_mass * radius))
        max_angvec = np.sqrt(self.__microthrust*np.pi/6. / (self.__sate_mass * (radius+8)))
        max_vt = max_angvec * radius

        delta_t = 10 * radius * np.pi/6 / (max_angvec * (radius+8))
        delta_vt = 2 * max_vt
        delta_vc = max_vt**2 * delta_t / radius
        delta_v = delta_vt + delta_vc
        #print(f"delta_vt: {delta_vt*1000:5.3f}, delta_vc: {delta_vc*1000:5.3f}")
        return delta_t/3600, delta_v


def main():
    darwin_sate = sate_motion(10e-3, 5e3)
    dis_list = [25, 50, 100, 225]
    for dis in dis_list:
        delta_t, delta_v = darwin_sate.linear_motion(dis)
        print(f"Distance: {dis:3d} m, delta_v: {delta_v*1000:5.3f} mm/s, delta_t: {delta_t:5.3f} hr")

    dis_list = [25, 50, 100, 250]
    for dis in dis_list:
        #max_angv, max_vt, delta_t = darwin_sate.circling_motion(dis)
        delta_t, delta_v = darwin_sate.circling_motion(dis)
        #print(f"Distance: {dis:3d} m, max_angv: {max_angv / np.pi * 180. * 3600:5.3f} arcsec/s, max_vt: {max_vt*1000:5.3f} mm/s")
        #print(f"Distance: {dis:3d} m, max_angv: {max_angv / np.pi * 180. * 3600:5.3f} arcsec/s, delta_t: {delta_t/3600:5.3f} hr")
        print(f"Distance: {dis:3d} m, delta_v: {delta_v*1000:5.3f} mm/s, delta_t: {delta_t:5.3f} hr")


if __name__ == "__main__":
    main()